Chapter 3: Imaging

---

3.6. Image Processing

Download

part of

MSE672: Introduction to Transmission Electron Microscopy

by Gerd Duscher, Spring 2022

Microscopy Facilities
Institute of Advanced Materials & Manufacturing
Materials Science & Engineering
The University of Tennessee, Knoxville

Background and methods to analysis and quantification of data acquired with transmission electron microscopes.

3.6.1. Load important packages

3.6.1.1. Check Installed Packages

import sys
from pkg_resources import get_distribution, DistributionNotFound

def test_package(package_name):
    """Test if package exists and returns version or -1"""
    try:
        version = get_distribution(package_name).version
    except (DistributionNotFound, ImportError) as err:
        version = '-1'
    return version

# pyTEMlib setup ------------------
if test_package('pyTEMlib') < '0.2021.3.22':
    print('installing pyTEMlib')
    !{sys.executable} -m pip install  --upgrade pyTEMlib -q
# ------------------------------
print('done')

done

3.6.1.2. Import all relevant libraries

# import matplotlib and numpy
#                       use "inline" instead of "notebook" for non-interactive plots
from PyQt5.QtGui import *
%pylab --no-import-all notebook
%gui qt
    
import warnings
warnings.filterwarnings('ignore')

# our blob detectors from the scipy image package
from skimage.feature import blob_dog, blob_log, blob_doh

# Multidimensional Image library
import scipy.ndimage as ndimage

import time
# Import libraries from the book
import sys
sys.path.insert(0, '../../pyTEMlib')
import pyTEMlib
import pyTEMlib.file_tools as ft          # File input/ output library
import pyTEMlib.kinematic_scattering as ks         # Kinematic sCattering Library



# it is a good idea to show the version numbers at this point for archiving reasons.
print('pyTEM version: ',pyTEMlib.__version__)

%pylab is deprecated, use %matplotlib inline and import the required libraries.
Populating the interactive namespace from numpy and matplotlib
You attempted to load a Gwyddion file, but this requires gwyfile.
     Please Load it with pip install gwyfile
You attempted to load a Nanonis file, but this requires Nanonispy.
     Please Load it with pip install nanonispy
You attempted to load a Nanonis file, but this requires Nanonispy.
     Please Load it with pip install nanonispy
Symmetry functions of spglib enabled
Using kinematic_scattering library version {_version_ }  by G.Duscher
pyTEM version:  0.2022.11.0

3.6.1.3. Load an atomic resolution image:

As an example we will use p1-3-hr3.dm3 in the example_data directory.

The image is of a sample of ZnO grown on graphene.

wid = ft.FileWidget()
wid

<pyTEMlib.file_tools.FileWidget at 0x7f79f24292a0>

# ---- Input ------
load_example = False
# -----------------
## Close file if one is open
try:
    dataset.h5_dataset.file.close()
except:
    pass

if load_example:
    file_name = '../example_data/p1-3-hr3.dm3'
else:
    file_name = wid.file_name 
print(file_name)    
# Load file

try:
    from PyQt5 import QtCore, QtWidgets, QtGui
except :
    print('Qt dialogs are not available')
    Qt_available = False

import sidpy
file_types = 'TEM files (*.dm3 *.emi *.ndata *.h5 *.hf5);;pyNSID files (*.hf5);;QF files ( *.qf3);;' \
                     'DM files (*.dm3 *.dm4);;Nion files (*.ndata *.h5);;All files (*)'
   
filename = sidpy.io.interface_utils.openfile_dialog_QT(file_types=file_types, file_path='.')

dataset = ft.open_file()
view = dataset.plot(vmax=0.1)

file type not handled yet.

---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
Input In [2], in <cell line: 2>()
      1 dataset = ft.open_file()
----> 2 view = dataset.plot(vmax=0.1)

AttributeError: 'NoneType' object has no attribute 'plot'

The X11 connection broke (error 1). Did the X11 server die?

3.6.2. Fourier Transform of Image

image = np.array(dataset-dataset.min())
fft_mag = np.abs((np.fft.fftshift(np.fft.fft2(image))))

## pixel_size in recipical space
rec_scale_x = 1/dataset.x[-1]  
rec_scale_y = 1/dataset.y[-1]

## Field of View (FOV) in recipical space please note: rec_FOV_x = 1/(scaleX*2)
rec_FOV_x = rec_scale_x * dataset.shape[0] /2.
rec_FOV_y = rec_scale_y * dataset.shape[1] /2.


## Field ofView (FOV) in recipical space
rec_extend = (-rec_FOV_x,rec_FOV_x,rec_FOV_y,-rec_FOV_y)

fig = plt.figure()
plt.imshow(np.log2(1+fft_mag).T, extent=rec_extend, origin = 'upper'); 
plt.xlabel('spatial frequency [1/nm]');

3.6.3. Shorter

from pyTEMlib import image_tools as it
fft_mag = it.power_spectrum(dataset)
fft_mag.plot()

3.6.4. More Contrast in Fourier Transform

Part of the low contrast is the high frequency noise in the diffractogram. Smoothing with a Gaussian (3pixels wide) increases the contrast on logarithmic intensity scale.

Additionally, the minimum and maximum intensity can be used to increase the plotted contrast. For that we want to exclude

• the center spot, which does not carry any spatial information and

• the high frequencies which have always low intensities.

# We need some smoothing (here with a Gaussian)
smoothing = 3
fft_mag2 = ndimage.gaussian_filter(fft_mag, sigma=(smoothing, smoothing), order=0)


#prepare mask for low and high frequencies
pixels = (np.linspace(0,image.shape[0]-1,image.shape[0])-image.shape[0]/2)* rec_scale_x
x,y = np.meshgrid(pixels,pixels);
mask = np.zeros(image.shape)

mask_spot = x**2+y**2 > 2**2 
mask = mask + mask_spot
mask_spot = x**2+y**2 < 10**2 
mask = mask + mask_spot

mask[np.where(mask==1)]=0 # just in case of overlapping disks

fft_mag3 = fft_mag2*mask
print(np.std(fft_mag2),np.mean(fft_mag2) ,np.max(fft_mag2) )
#print(np.std(fft_mag2[np.where(mask==2)]),np.mean(fft_mag2[np.where(mask==2)]) ,np.max(fft_mag2[np.where(mask==2)]))


minimum_intensity = np.log2(.1+fft_mag2)[np.where(mask==2)].min()*0.95
#minimum_intensity = np.mean(fft_mag3)-np.std(fft_mag3)
maximum_intensity = np.log2(.1+fft_mag2)[np.where(mask==2)].max()*1.05
#maximum_intensity =  np.mean(fft_mag3)+np.std(fft_mag3)*2

print(minimum_intensity,maximum_intensity)
fig = plt.figure()
plt.imshow(np.log2(.1+fft_mag2).T, extent=rec_extend, origin = 'upper',vmin=minimum_intensity, vmax=maximum_intensity)
plt.xlabel('spatial frequency [1/nm]');

0.5914227754441985 11.579143555569678 18.297265113724112
3.340085191885885 4.176874317680241

3.6.5. Spot Detection

This diffractogram is now good enough to detect the spots with a blob detector.

# ---- Input ---------
spot_threshold = .4
# --------------------

## Needed for conversion from pixel to Reciprocal space
rec_scale = np.array([rec_scale_x, rec_scale_y,1])
center = np.array([int(image.shape[0]/2), int(image.shape[1]/2),1] )

fft_mag2 = ndimage.gaussian_filter(fft_mag, sigma=(smoothing, smoothing), order=0)
fft_mag2 = fft_mag2-fft_mag2.min()
#fft_mag2 = fft_mag2/fft_mag2.max()

# spot detection ( for future reference there is no symmetry assumed here)

spots_random =  (blob_log(fft_mag,  max_sigma= 5 , threshold=spot_threshold)-center)*rec_scale

print(f'Found {spots_random.shape[0]} reflections')
spots_random[:,2] = np.linalg.norm(spots_random[:,0:2], axis=1)
spots_index = np.argsort(spots_random[:,2])
spots = spots_random[spots_index]

## plot Fourier transform and found spots
fft_mag.plot()
plt.gca().scatter(spots[:,0], spots[:,1], c='Red',  alpha = 0.2, label='spots');
#plt.xlim(-0.2,0.2)
#plt.ylim(-0.2,0.2);

Found 23 reflections

3.6.6. Adaptive Fourier Filtering

We mask the fourier transformed image so that the information can pass through is selected.

The information is in the spots and in the center of the Fourier transformed image, the rest is noise.

Please modify the radius of the mask of the reflections and the low-pass area in the code below and notice the effects on the Fourier filtered image.

# Input 
reflection_radius = 0.3 # in 1/nm
low_pass = 1/1.34 # in 1/nm diameter of mask for low pass filter
FOV_x = dataset.x[-1]
FOV_y = dataset.y[-1]


#prepare mask
pixels = (np.linspace(0,image.shape[0]-1,image.shape[0])-image.shape[0]/2)* rec_scale_x

x,y = np.meshgrid(pixels,pixels);
mask = np.zeros(image.shape)

# mask reflections
for spot in spots:
    mask_spot = (x-spot[0])**2+(y-spot[1])**2 < reflection_radius**2 # make a spot 
    mask = mask + mask_spot# add spot to mask
    
# mask zero region larger (low-pass filter = intensity variations)

mask_spot = x**2+y**2 < low_pass**2 
mask = mask + mask_spot

mask[np.where(mask>1)]=1 # just in case of overlapping disks

plt.figure()
ax1 = plt.subplot(1,2,1)
#ax1.imshow(mask)
fft_filtered = np.fft.fftshift(np.fft.fft2(image))*mask.T
ax1.imshow(np.log(1+np.abs(fft_filtered)).real,extent=rec_extend, origin = 'upper')
plt.xlabel('reciprocal distance [1/nm]')
ax2 = plt.subplot(1,2,2)
filtered = np.fft.ifft2(np.fft.fftshift(fft_filtered))

real_extent = (0,FOV_x,FOV_y,0)
    
ax2.imshow(filtered.real,extent=real_extent, origin = 'upper')
plt.xlabel('distance [nm]');

3.6.6.1. Plot the image a bit bigger

plt.figure()
plt.imshow(filtered.real, extent=real_extent, origin = 'upper')
plt.xlabel('distance [nm]');

3.6.7. Check on filtered images

We don’t want to filter anything out that caries information, or at least we want to be aware of that. An easy check is to subtract the filtered imag fromt he image and to determine that only noise is left.

Please note that any processing can be easily determined in the Fourier transformed, so be meticulous on reporting what you did to an image.

filtered = np.array(filtered)
plt.figure()
ax1 = plt.subplot(1,2,1)
im1 = ax1.imshow(image-filtered.real,extent=real_extent, origin = 'upper', vmin=-100, vmax=100)
plt.xlabel(' distance [nm]')
plt.colorbar(im1)
ax2 = plt.subplot(1,2,2)
fft_difference = np.fft.fftshift(np.fft.fft2(image-filtered.real).T)
im2 = ax2.imshow(np.log(1+np.abs(fft_difference)).real.T,extent=rec_extend, origin = 'upper')
plt.xlabel('reciprocal distance [1/nm]');

3.6.8. Conclusion:

Fourier Transform is used to evaluate image quality and Fourier Filtering is used to eliminate high frequency noise.

If the computer is used to find the spots in the diffractogram, we do not have to select them ourselfs. A carefull selection will,however, yield the same result.

3.6.9. Navigation

• Up Chapter 3: Imaging

• Back: Defocus Thickness Map

• Next: Image Analysis

• List of Content: Front